% 数值雅可比矩阵计算

function J = Jacobian(f, x, eps)
    n = length(x);
    fx = f(x);
    m = length(fx);
    J = zeros(m,n);
    for i = 1:n
        dx = zeros(n,1);
        dx(i) = eps;
        J(:,i) = (f(x + dx) - fx) / eps;
    end
end